Exact time-dependent density functional theory for impurity models 
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We employ the density matrix renormalization group to construct the exact time-dependent ex- 
change correlation potential for an impurity model with an applied transport voltage. Even for 
short-ranged interaction we find an infinitely long-ranged exchange correlation potential which is 
built up instantly after switching on the voltage. Our result demonstrates the fundamental diffi- 
culties of transport calculations based on time-dependent density functional theory. While formally 
the approach works, important information can be missing in the ground-state functionals and may 
be hidden in the usually unknown non-equilibrium functionals. 



I. INTRODUCTION 

In recent years a combination of Kohn-Sham density 
functional theory (DFT) and the Landauer approach 
to transport has been developed that allows an ab- 
initio calculation of the current- voltage characteristics of 
molecules that are attached to reservoirs. Early com- 
parisons between calculated and experimental conduc- 
tances yielded discrepancies of several orders of magni- 
tude, while in more recent comparisons a typical discrep- 
ancy of one order of magnitude was reported!^ Concern- 
ing experiments, the reproducibility of I-V characteris- 
tics in molecular electronics poses a problem, as the con- 
tact configuration may change from sample to sample, 
which makes the comparison between theory and exper- 
iment difficult. 

First DFT studies of lattice models like the Hubbard 
model or spinless fermions date already back to the late 
1980sPH At that time ground state DFT was the main 
focus of interest. In the 2000s a lattice version of the 
local density approximation was constructed using the 
Bethe ansatz solution of the Hubbard model, and the ac- 
curacy of the method was tested in some detailPH Due 
to the mentioned discrepancies between measured and 
calculated conductances through molecules, recently the 
focus turned to the DFT description of transport. In the 
conventional DFT approach to the two-terminal linear 
conductance of a molecular system the Kohn-Sham equa- 
tions are solved in order to obtain the electronic struc- 
ture, and the conductance of the Kohn-Sham system, 
Gks j is used as a theoretical estimate of the true conduc- 
tance, G. There are thus two possible sources for errors, 
(i) the electronic structure calculation and (ii) the identi- 
fication of G with Gks ■ While in complex realistic struc- 
tures it is not possible to quantify (i) and (ii) this is possi- 
ble for lattice models, as demonstrated by Schmitteckert 
and EversP They observed that the Kohn-Sham conduc- 
tance, GkS: based on the exact ground state exchange 
correlation potentials becomes accurate, i.e. Gks = G, 
close to (isolated) conductance resonances. Remarkably 
this holds even if the spectral properties of the Kohn- 
Sham system strongly differ from the true onesP^ For 
example, the zero temperature conductance through a 



Kondo impurity is captured correctly by Gks- The rea- 
son for this coincidence is the Friedel sum rule which 
guarantees that G and Gks hav e identi cal functional de- 
pendencies on the charge density. ^! 11 ! 12 ! However, as soon 
as the conductance is not given by a local Friedel sum 
rule there can be ord ers of m agnitude between Gks an d 
the true conductance) 10 * 13 -^ and it can even be paramet- 
rically wrong.^ It is important to note that this is not 
a general failure of DFT; rather the assumption that the 
conductance of the physical electrons is given by the con- 
ductance of the Kohn-Sham particles breaks down. It has 
been stated several times in the literature that this dis- 
crepancy can be corrected by including dynamic contri- 
butions to the exchange-correlation potential that are not 
captured by ground state functionals pHES p or th e linear 
current / through an impurity this dynamic exchange- 
correlation potential V xc renormalizes the voltage, such 
that 

I = GV = G KS (V + V* C )- (1) 

It is not clear whether a similar voltage renormalization 
is also present in the non-linear response. 

The goal of our study is to investigate the nature of 
the dynamical response within a DFT framework beyond 
the linear transport regime. To this end we construct 
the exact time-dependent exchange-correlation potential 
for an impurity model with an applied transport voltage. 
Specifically we construct the time-dependent Kohn-Sham 
potentials for a two-terminal transport setup, which by 
construction yields the correct physical current within 
the time-dependent DFT description. In this paper 
we will concentrate on a one-dimensional lattice model, 
where we can directly compare time-dependent DFT with 
accurate numerical and analytical results obtained by 
many-body techniquesP In order to obtain meaningful 
results we have to solve much larger systems (here: 240 
lattice sites) than for instance in Ref. [TS] (6 to 12 lat- 
tice sites). In addition the densities have to be calcu- 
lated with an accuracy better than ~ 10~ 6 ; otherwise 
the reverse engineering of the time-dependent exchange- 
correlation potentials fails. 

We consider the interacting resonant level model, i.e. 
a one-dimensional model of spinless fermions, where a 
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single interacting level is coupled to a left and a right 
lead. 
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Here t = 1 is the hopping amplitude, Nr and ATr, are the 
numbers of sites in the left and right lead, and U is the 
interaction on the contact link. All the data presented in 
this article have been obtained for the half-filled lattice 
model. We assume that at time T < the system is 
in the ground state. At T = we include a voltage 
drop by applying a potential eV/2 (— eV/2) on the left 
(right) lead with a linear crossover on a scale of several 
sites left and right of the impurity. Then we follow the 
time evolution of the system. On time scales that are 
shorter than the transit time T t = L\ eg _&/vp (vp — 2t is 
the Fermi velocity) the finite leads act as reservoirs, such 
that a time-dependent simulation allows to extract the 
current-voltage relation corresponding to infinite leads, 
see Refs. H] and [20] for details. Note that in our m odel 
the linear conductance is given by a Friedel sum rulel^il 
Therefore, V xc vanishes in the linear regime. 

We consider a Hamiltonian H with parameters U and 
t' and a Kohn-Sham Hamiltonian Hks with the same 
structure as Eq. ^ using different parameters Uks = 
and t' KS . The essence of DFT is a one-to-one correspon- 
dence between local densities and potentials: Starting 
with an uniquely defined initial state there exists an - 
up to a gauge transformation - unique set of potentials 
{vks(T)} such that the time-dependent densities {n(T)} 
are identical for H and 7/ks- The numerical task is to 
calculate first the time-dependent local density for H and 
in the second step the set of potentials {«ks} f° r -Hrs- 
We fix the gauge by imposing that the sum of the poten- 
tials is zero. To calculate the particle density we follow 
two strategies: (a) The particle density in an interacting 
system is obtained through time-dependent density ma- 
trix renormalization group (td-DMRG) as described in 
Refs. HUJ and Since these calculations are extremely 
time-consuming we also (b) study a toy model of non- 
interacting fermions, where the time-dependent density 
can be obtained straightforwardly by exact diagonaliza- 
tion. As mentioned before, when calculating {i>ks} from 
the densities a high accuracy is necessary. We use an it- 
erative procedure that stops when the densities in H and 
in JTks agree within an error of 10 -10 . 




FIG. 1: The Kohn-Sham potential vks as a function of time 
for the first ten lattice sites from the left edge of the chain. 
The two data sets are for a total chain length N = 200 and 
N = 400. For comparison the potential v — 0.2t in the toy 
model (corresponding to a voltage eV = OAt) is also shown 
as dotted line. 



II. NONINTERACTING TOY MODEL 

We start by considering the non-interacting version of 
the Hamiltonian ([2]), i.e. we set the interaction U as well 
as Uks to zero. There remains only one free parame- 
ter, the hopping amplitude t' between the impurity site 
and the leads. In our toy model we choose t' — 0.5t 
and t' KS = 0.3i. Typical results for the time-dependent 
Kohn-Sham potentials are shown in Fig. [T] In the toy 
model a potential eV = OAt corresponding to a potential 
v = ±0.2i in the left (right) lead is switched on at T = 0. 
As a response a current flows through the resonant level 
that becomes stationary after T ~ 30 (the time is in 
units h/i). By definition, the same current has to flow 
also in the Kohn-Sham system; here we find initially a 
time-dependent voltage which becomes stationary after 
T « 30. Remarkably the exchange-correlation potential 
(the difference between «ks and v) is nearly position in- 
dependent in the leads, and appears immediately after 
switching on the voltage. In the figure, data for a chain 
of 200 lattice sites and a chain of 400 lattice sites are 
presented. In both cases we plot the potentials for the 
first ten sites at the left edge of the chain. For the longer 
chain vks remains stationary until the end of the sim- 
ulation, and only a single line is seen meaning that the 
potential is homogeneous in space. In the shorter chain 
the potential becomes position and time-dependent after 
T sa 50. This happens since the simulation time is longer 
than the transit time T t = 100/2 = 50. 

The Kohn-Sham potential as a function of position is 
illustrated in Fig. [2] demonstrating that a major effect of 
the exchange-correlation potential is a time-independent 
renormalization of the local potentials, which implies an 
additional voltage as anticipated in Eq. ([IJ. A reason- 
able estimate for the voltage renormalization is obtained 
by comparing the I-V characteristics for two noninter- 
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FIG. 2: Local potential Vi in the toy model and in the Kohn- 
Sham system close to the impurity (i = 0) at time T = 12 
(solid curve) and T = 20 (symbols). 
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FIG. 3: Current (in units of et/h) as a function of voltage 
for two non-interacting systems with t' = 0.3t and t' = 0.5i, 
respectively. In order to generate the same current in both 
systems (e.g. I = 0.4 as indicated by the dotted line) different 
potentials, eVo.3 and eVb.s, have to be applied. The voltage 
renormalization V0.3/V0.5 displayed in the inset gives a rea- 
sonable estimate of the voltage renormalization observed in 
the Kohn-Sham Hamiltonian for the toy model. 



acting models with t' — 0.3t and t' — 0.5t, see Fig. [3] 
The current I = 0.3883 (in units of et/h) corresponds to 
eV = OAt for t' = 0M and eV = 0.6061i for t' = 0M 
which is close to the voltage we found in the Kohn-Sham 
system, see Fig.[l] Notice that the maximum current that 
can be achieved is larger in the case t' — 0.5t than for 
t' = 0.3t. For DFT this implies that no stationary Kohn- 
Sham potential can generate such large currents. The- 
oretically, time-dependent Kohn-Sham potentials could 
generate a stationary current, however, the proofs of ex- 
istence of time-dependent Kohn-Sham potential^HSS] q\ 
not apply to the present situation, so that possibly the 
high voltages are not w-representable, compare Refs. |26| 
and [27l A similar situation will also arise in the next 
section where we investigate the model with interaction. 



III. DMRG RESULTS FOR THE INTERACTING 
MODEL 

We now turn our attention to the interacting case 
and choose H with U — 2t and t' = 0.3t, whereas 
for the Kohn-Sham Hamiltonian we set Uks = and 
t' KS = 0.3t. The I-V characteristics of this model is 
known analytically^ from field theoretical methods. For 
example, the current is given in closed form by^ 
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where 3F2 is a hypergeometric function and eV c — rt' 4 ^ 3 
sets the scale separating a charge 2e dominated low volt- 
age regime, and a charge e/2 dominated high voltage 
regime.^ The regularizati on r ~ 3.2 was determined 
from numerical simulation! 2 * 28 * 29 ^ The study of the toy 
model in the previous section shows that in order to ob- 
tain meaningful results the simulation time should be of 
the order T « 30 or longer, so that the length of each 
lead should be at least L lcad rj 60 sites. This means 
that the excitations originating form the charge quench 
at the impurity should not reach the boundary of the 
leads within the simulation time. In order to obtain the 
time-dependent densities we apply the full td-DMRG 22 
as it has already proven to provide accurate results that 
agree perfectly with the analytic solution. 

The quench with a symmetric charge imbalance leads 
to oscillations with frequency eV/2 during the transient 
time and finite size induced Josephson like oscillations 
with frequency eV in the steady state regime PS There- 
fore, when the effective voltage of the Kohn-Sham system 
is different from the voltage of the physical one, the time- 
dependent potentials have to compensate this mismatch. 
In order to reduce these finite time and finite size effects 
we ramp up the voltage linearly from T = Q . . . 5, leading 
to reduced amplitudes of the transient oscillations. Fig. [4] 
shows the Kohn-Sham potential in the left lead. As in 
the non-interacting toy model the main effect is a voltage 
renormalization. 

To obtain numerically converged data is a hard task. 
The figure shows results from three different DMRG runs 
varying the numerical accuracy. In all three runs the 
time-dependent particle density is almost identical, varia- 
tions occur on the scale 10 -6 . The Kohn-Sham potentials 
turn out to be very sensitive functions of the particle den- 
sity and show pronounced differences in the three cases. 
Only for a very large number of states per block in the 
DMRG calculations (highest accuracy) the final result 
becomes smooth. In Fig. [5] time-dependent Kohn-Sham 
potentials for different voltages are depicted. Whereas for 
a small voltage we are able to find Kohn-Sham potentials 
for all times, this is not the case for large voltage where 
the Kohn-Sham potential diverges at finite T. It will be 
interesting to see in future research, whether the time 
scale of the singularity has a deeper meaning. Since the 
w-representability is not guaranteed in the lattice model, 
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FIG. 4: Local potentials on the first ten sites of the left lead 
as a function of time for an interacting system (U = 2t, 240 
lattice sites). Here we smoothly switched on the external volt- 
age V = OAt between T = and T = 5. The number of states 
per block kept in the DMRG calculations varied between 2000 
and 4000. 
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FIG. 5: Time-dependent Kohn-Sham potentials on the first 
ten sites of the left lead varying the voltage applied to the 
system (t' = 0.3i and U — 2t). While for a small voltage 
(eV = OAt) we are able to find a set of Kohn-Sham potentials 
for all times, for larger voltages (eV = 0.6t, eV = 0.8t) this 
is only possible for short times. The inset shows the I-V 
characteristics for U = and U = 2t. The curve for U = 2t 
is given by Eq. (Ky and the symbols are the DMRG results. 



we believe that Kohn-Sham potentials do not exist for all 
times in these cases, compare also the discussion relating 
to Fig. U 

In the inset of Fig. [5] we show the I-V characteristics 
of the noninteracting (U — 0) and interacting (U = 2t) 
resonant level model. Again, there is a regime, where the 
current in the interacting case is higher than any current 
achievable in the noninteracting case. Since the noninter- 
acting case is similar to the Kohn-Sham system, this is a 
hint that no stationary Kohn-Sham potential might exist 
in this regime. The cases where we find diverging Kohn- 
Sham potentials are indeed in the regime where the cur- 
rent for U — 2t exceeds the one for U = 0. Of course, the 
Kohn-Sham systems at finite voltage have additional po- 
tentials with a spatial structure which is not captured by 
a single number, V KC . It is also important to note, that 
the Kohn-Sham potentials should not create a current 
outside the "light cone" vfT. The continuity equation 
then implies that any oscillations of the potentials out- 
side this regime have to be instantaneous and constant in 
space including the reservoirs. While one can gauge the 
resulting V xc into a time-dependent phase exp(iF xc T) of 
a hopping element at the impurity, one loses the prop- 
erty of an at least locally stationary system found in this 
work. 



IV. SUMMARY 

We calculated the exchange-correlation potential for 
a one-dimensional model system with an applied trans- 
port voltage. Specifically we considered an impurity 
model with short-ranged interaction and non-interacting 
leads. Immediately after switching on the voltage the 
exchange-correlation potential appears deep inside the 
leads. This long-ranged potential is purely dynamic, and 
therefore approximations based on equilibrium function- 
als and short-ranged approximations have to fail. 
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